• Todennäköisyys ja Bayesilainen tilastotiede
  • Bayesin kaava
  • Uskottavuusfunktio
  • Priorijakauma
  • Posteriorijakauma
• Bayesilaiset menetelmät käytännössä
• Kalvot löytyvät osoitteesta github.com/hpesonen

• Bayesiläisen tilastotieteen keskiössä on maailman mallittaminen
• Bayesilainen tilastotiede mallittaa todennäköisyyden epävarmuutena todellisuudesta
• 0 ja 1 ovat skaalan kaksi ääripäätä ja vastaavat täyttä varmuutta tapahtumasta
  • Epävarmuus on subjektiivinen käsite
  • Maailman havannointi lisää informaatiota todellisuudesta ja vähentää epävarmuutta
• Bayesin teoreeman käyttö ei tee kenestäkään Bayesilaista, mutta epävarmuuden ilmaiseminen todennäköisyyden avulla tekee

• 0 ja 1 ovat skaalan kaksi ääripäätä ja vastaavat täyttä varmuutta tapahtumasta

• Informaation lisäämisen mekanismi ilmaistaan Bayesin teoreeman avulla
  • Bayesin teoreema on siis työkalu, jota käytetään Bayesilaisessä analyysissä

\[p(\theta \mid y) = \frac{p(y \mid \theta)p(\theta)}{p(y)}\]

\[p(\theta \mid y) = \frac{p(y \mid \theta)p(\theta)}{p(y)}\]

• \(\theta\): Mallin parametrit
• \(y\): Data
• \(p(\theta \mid y)\): Posteriorijakauma
• \(p(y \mid \theta)\): Uskottavuusfunktio
• \(p(\theta)\): Priorijakauma
• \(p(y)\): Evidenssi/Marginaali

• Kolikonheitto on klassinen esimerkki informaation kertymisestä havaintoaineiston kasvaessa
• Tarkastellaan kolikkoa ja halutaan selvittää onko kolikko reilu
  • Merkitään parametrilla \(p\) todennäköisyyttä saada kolikonheitossa klaava
• Voidaan ajatella, että kolikonheitto on Bernoullin koe
  • Klaava (\(y=1\)) saadaan todennäköisyydellä \(\theta\), jolloin kruunan todennäköisyys on \(1-\theta\)
  • Uskottavuusfunktio Bernoullin kokeelle on

\[p(y \mid \theta) = \theta^y(1-\theta)^{1-y}\]

• Toistettaessa riippumattomasti kolikonheitto-koetta \(N\) kertaa uskottavuusfunktio voidaan ilmaista muodossa

\[p(y_1, \ldots, y_N \mid \theta) = \prod_{i=1}^N p(y_i \mid \theta) = \theta^{\sum_{i} y_i} (1- \theta)^{N-\sum_{i}y_i} \]

• Miten mallintaisimme epävarmuuttamme ennen ensimmäistäkään havaintopistettä?
  • Minkälainen priorijakauma vastaisi epävarmuuttamme klaavan esiintymisestä?

• Uskoaksemme on hyvin epätodennäköistä, että melko varmasti saataisiin aina \(0/1\)
• Uskoaksemme kolikko on melko todennäköisesti reilu (\(\theta=0.5\)), mutta emme voi poissulkea mahdollisuutta, että se ei ole
• Mallinnamme prioritietämystä jakaumalla \(p(\theta) \propto \theta^4(1-\theta)^4\) (\(\mathsf{Beta}(5,5)\)-jakauma)

• Tarkastellaan miten posteriorijakauma muuttuu heittosarjan koon kasvaessa
• Merkitään dataa \(y = (n,N-n)\) :
• Esim. \(y = (10,5)\) vastaa \(15\) heiton sarjaa, joista \(10\) klaavaa ja \(5\) kruunaa
• Tässä tapauksessa posteriori-jakauma voidaan ilmaista analyyttisesti

\[ \begin{aligned} p(\theta \mid y = (n,N-n)) & \propto p(y=(n,N-n) \mid \theta) p(\theta) \\ & = \theta^{n}(1-\theta)^{N-n}\theta^4(1-\theta)^4 \\ & = \theta^{n+4}(1-\theta)^{N-n+4} \end{aligned} \]

• Havaintojen määrän kasvaessa posteriori-jakauma lähestyy uskottavuusfunktiota
  • Priori-todennäköisyyden merkitys vähenee

• Tulkitaan posteriorijakauma

• Frekventistisen tilastotieteen luottamusväliä vastaava käsite Bayesilaisessa tilastotieteessa on todennäköisyysväli
• Määritelmän mukaisesti \(\alpha\)-todennäkäisyysväli (yleisemmin todennäköisyysalue) on mikä tahansa parametrijoukko, joka sisältää \(\alpha\) todennäköisyyttä \[p(\theta \in \Theta \mid x) = \alpha\]
• Käytännössä kyse on integraalista (tai summasta diskreetin satunnaismuuttujan tapauksessa) \[p(\theta \in \Theta \mid x) = \int_\Theta p(\theta \mid x) ~\mathsf{d} \theta = \alpha\]

• Todennäköisyysvälit eivät ole yksikäsitteisiä
  • Pienin väli/alue, joka sisältää \(\alpha\) todennäköisyyden
  • Väli, jonka alapuolelle jäävä todennäköisyys vastaa yläpuolelle jäävä todennäköisyyttä
  • Odotusarvon ollessa olemassa, väli joka on keskitetty odotusarvoon

• \(0.95\)-Todennäköisyysväli lyhenee havaintojen määrän kasvaessa

• Ero frekventistiseen luottamusväliin on, että Bayesilaisen tulkinnan mukaisesti parametrin arvo kuuluu \(\alpha\)-todennäköisyysväliin todennäköisyydellä \(\alpha\)

• Mallien vertailu tapahtuu mallien posterioritodennäköisyyden avulla
• Jotta mallin \(M_i\) todennäköisyys voidaan laskea, tarvitaan malleille prioritodennäköisyydet \(p(M_i)\) kuten mille tahansa parametrille

\[\sum_{i=1}^m p(M_i) = 1\]

• Mallien posterioritodennäköisyydet ovat

\[\sum_{i=1}p(M_i \mid y) = 1\]

• Oletetaan, että mallin \(M_i\) määrittää sen parametrit \(\theta_i\)
• Tällöin mallin posterioritodennäköisyys lasketaan \[ p(M_i \mid y) = \frac{p(y \mid M_i) p(M_i)}{p(y)} \]
• Nyt mallin uskottavuus, evidenssi, lasketaan \[p(y \mid M_i) = \int p(y \mid \theta_i, M_i) ~ p(\theta_i \mid M_i) ~ \mathsf{d} \theta_i\]

• Nykyisin käytetään usein mallien vertailuun posteriorijakauman laskemisen jälkeen Leave-one-out ristiinvalidointia
• Laskettava suure on \[\text{elpd}_\text{loo} = \sum_{i=1}^N \log p(y_i \mid y_{-i}) = \sum_{i=1}^N \log \int p(y_i \mid \theta)p(\theta \mid y_{-i}) \mathsf{d}\theta\]

• Suureen \(\text{elpd}_\text{loo}\) laskemissa käytetään hyväksi prediktiivisiä jakaumia
• Prediktiivisiä jakaumia voidaan yleisemminkin hyödyntää havaintomallien sekä priorijakaumien järkevyyden tarkasteluun
• Nämä voidaan jakaa prioriprediktiivisiin sekä posterioriprediktiivisiin jakaumiin

• Kun havaintomalli sekä tämän parametrien priorijakauma on mallitettu, minkälaista havaintoaineistoa malli ennustaa

\[p(y) = \int p(y \mid \theta) p(\theta) \mathsf{d} \theta\]

• Usein prioritieto liittyy nimenomaisesti tietämykseen mallin tuottamista realistisista arvoista

• Prioriprediktiivinen malli selventää parametrien priorijakauman vaikutusta havantoihin
  • Tuottavatko jotkut parametrikombinaatiot mahdottomia/hyvin epätodennäköisiä havaintoja?

• Kun havaintomalli sekä tämän parametrien priorijakauma on mallitettu ja posteriorijakauma on laskettu havaintoaineistoa käyttäen, minkälaista uutta havaintoaineistoa \(y^\prime\) malli ennustaa

\[ \begin{aligned} p(y^\prime \mid y) & = \int p(y^\prime \mid \theta, y) p(\theta \mid y) \mathsf{d} \theta \\ & = \int p(y^\prime \mid \theta) p(\theta \mid y) \mathsf{d} \theta \end{aligned} \]

• Posterioriprediktiivistä mallinnusta voidaan käyttää tarkastelemaan opittujen mallien järkevyyttä

• Posteriori-jakauman laskeminen analyyttisessä muodossa onnistuu ainoastaan erikoistapauksissa
• Käytännössä laskeminen tapahtuu muodostamalla otosparvi posteriorijakaumasta
  • Importance sampling
  • Sequential Monte Carlo
  • Markov Chain Monte Carlo
  • Hamiltonian Monte Carlo

• Otosparven avulla voidaan likimääräisesti laskea \(\mathsf{E}(g( \theta)\mid y)\) mille tahansa funktiolle \(g(\cdot)\)

• Nykyisin otoksen muodostamiseen on yksinkertaista laajaltikäytettyjen ohjelmistojen ansiosta
• Suosituin ohjelmistoista on Stan
  • Stan
    • PyStan
    • RStan
• rstanarm : Bayesian applied regression modeling
  • R:stä tuttu käyttöliittymä Bayesilaisesta analyysista kiinnostuneille

• Data : kidiq
• Muuttujat
  • Lapsen ÄO
  • Äidillä high school-tutkinto indikaattorimuuttuja
  • Äidin ÄO
  • Äidin ikä

## 
## Model comparison: 
## (ordered by highest ELPD)
## 
##       elpd_diff se_diff
## post4   0.0       0.0  
## post3  -3.5       2.5  
## post2  -6.0       3.9  
## post1 -42.2       8.7

## 
## Model comparison: 
## (negative 'elpd_diff' favors 1st model, positive favors 2nd) 
## 
## elpd_diff        se 
##      42.2       8.7

## 
## Model comparison: 
## (negative 'elpd_diff' favors 1st model, positive favors 2nd) 
## 
## elpd_diff        se 
##       3.5       2.5

## 
## Model comparison: 
## (negative 'elpd_diff' favors 1st model, positive favors 2nd) 
## 
## elpd_diff        se 
##       6.0       3.9

## 
## Call:
## glm(formula = lot1 ~ log(u), family = Gamma, data = clotting)
## 
## Deviance Residuals: 
##      Min        1Q    Median        3Q       Max  
## -0.04008  -0.03756  -0.02637   0.02905   0.08641  
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -0.0165544  0.0009275  -17.85 4.28e-07 ***
## log(u)       0.0153431  0.0004150   36.98 2.75e-09 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for Gamma family taken to be 0.002446059)
## 
##     Null deviance: 3.51283  on 8  degrees of freedom
## Residual deviance: 0.01673  on 7  degrees of freedom
## AIC: 37.99
## 
## Number of Fisher Scoring iterations: 3

## 
## Call:
## glm(formula = lot2 ~ log(u), family = Gamma, data = clotting)
## 
## Deviance Residuals: 
##      Min        1Q    Median        3Q       Max  
## -0.05574  -0.02925   0.01030   0.01714   0.06371  
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -0.0239085  0.0013265  -18.02 4.00e-07 ***
## log(u)       0.0235992  0.0005768   40.91 1.36e-09 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for Gamma family taken to be 0.001813354)
## 
##     Null deviance: 3.118557  on 8  degrees of freedom
## Residual deviance: 0.012672  on 7  degrees of freedom
## AIC: 27.032
## 
## Number of Fisher Scoring iterations: 3

## stan_glm
##  family:       Gamma [inverse]
##  formula:      clot_time ~ log_plasma * lot_id
##  observations: 18
##  predictors:   4
## ------
##                    Median MAD_SD
## (Intercept)        -0.016  0.007
## log_plasma          0.015  0.003
## lot_id2            -0.007  0.013
## log_plasma:lot_id2  0.008  0.006
## 
## Auxiliary parameter(s):
##       Median MAD_SD
## shape 7.942  2.577 
## 
## Sample avg. posterior predictive distribution of y:
##          Median MAD_SD
## mean_PPD 31.944  4.859
## 
## ------
## * For help interpreting the printed output see ?print.stanreg
## * For info on the priors used see ?prior_summary.stanreg

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

## (Intercept)     dist100     arsenic 
##       0.000      -0.894       0.461

## Loading required package: gridExtra

## 
## Attaching package: 'gridExtra'

## The following object is masked from 'package:dplyr':
## 
##     combine

## Warning: Removed 39 rows containing missing values (geom_point).